
#----------------------------------------------------------------------------------------------------------------#
#
# This code replicates "Figure 1: Average real monthly expenditures per EGM" that appears 
# in the paper "How windfall income increases spending at poker machines".
# 
# This code requires a shape file from the Australian Bureau of Statistics. 
#
# Created 21 Dec.c 2012 by Kyle Peyton
# Last update 17 Jan. 2012 by Kyle Peyton
#
# Questions, comments, complaints? contact Kyle Peyton: kyle.peyton@unimelb.edu.au 
#----------------------------------------------------------------------------------------------------------------#


rm(list=ls())

#packages required to replicate this figure. You will need to download these if not already installed. 
require('colorspace')
require("ggplot2")  
require("rgdal") 
require("maptools")  
require('gpclib')
require('foreign')
gpclibPermit() 


#----------------------------------------------------------------------------------------------------------------#
#read in shape file from the ABS (Australian Burea of Statistics) website  
#   shape files can be obtained here: 
#       http://www.abs.gov.au/AUSSTATS/abs@.nsf/DetailsPage/1259.0.30.001July%202011?OpenDocument 
#----------------------------------------------------------------------------------------------------------------#

#import shape file from ABS; path name depends on your system. 
shp <- readShapeSpatial( fn = '/Users/Kyle/Dropbox/Pokies/Replication/LGA11aAust', verbose = TRUE, delete_null_obj=TRUE)  
summary(shp)  

#list of names in the shape file
list_names = names(shp)  

# create a subset of the spatial data, take only the LGAs for the state of Victoria  (e.g. STATE_CODE ==2)
vic <- shp[which(shp$STATE_CODE==2),]  

#fortify data so ggplot can map it.  
u <- unique(vic$LGA_CODE11)        
out <- fortify(shp[u,], region="LGA_CODE11")  

# read in Pokies data compiled from Victorian Commission for Gambling and Liquor Regulation (VCGLR)
pokies <- read.dta("/Users/Kyle/Dropbox/Pokies/Replication/pokies01.dta") 

#subset real monthly expenditure per EGM along with key variables
pokiesubset <- subset(pokies, select=c(LGA_short, LGA_ABScode, Month, Year, rMonExp_EGM))

#create datafame of average expenditures by LGA over July 2004 - June 2012
pokiesavg <- aggregate(pokiesubset$rMonExp_EGM, by=list(LGA_short=pokiesubset$LGA_short,
                                                        LGA_ABScode=pokiesubset$LGA_ABScode), 
                       FUN=mean)

colnames(pokiesavg)[3] <- 'rMonExp_EGM' 

#eliminate '*' on the LGA codes in preparation for merge. 
pokiesavg$LGA_ABScode <- sub("\\*","", pokiesavg$LGA_ABScode)

#Increase memory to 3GB
memory.limit(size=300000)

# merge pokie expenditures with spatial data  
m <- merge(vic, pokiesavg,  by.x="LGA_CODE11", by.y="LGA_ABScode", all.x=TRUE)  

#Force subsets to aggregate LGA codes. Note that the VCGLR aggregates data for
#these LGAs. This must be done to match their estimates. 

#pyrenees to Ararat
m$rMonExp_EGM[m$LGA_CODE11 == 25990] <- m$rMonExp_EGM[m$LGA_CODE11 == 20260]

#Loddon (23940), Mount Alexander (25430) and Hepburn (22910) to Central Goldfields
m$rMonExp_EGM[m$LGA_CODE11 == 23940] <- m$rMonExp_EGM[m$LGA_CODE11 == 21670]
m$rMonExp_EGM[m$LGA_CODE11 == 25430] <- m$rMonExp_EGM[m$LGA_CODE11 == 21670]
m$rMonExp_EGM[m$LGA_CODE11 == 22910] <- m$rMonExp_EGM[m$LGA_CODE11 == 21670]

# Moyne (25490), Southern Grampians (26260) and Queenscilffe (26080) to Corangamite (21830)
m$rMonExp_EGM[m$LGA_CODE11 == 25490] <- m$rMonExp_EGM[m$LGA_CODE11 == 21830]
m$rMonExp_EGM[m$LGA_CODE11 == 26260] <- m$rMonExp_EGM[m$LGA_CODE11 == 21830]
m$rMonExp_EGM[m$LGA_CODE11 == 26080] <- m$rMonExp_EGM[m$LGA_CODE11 == 21830]

#Towong(26670) to Mansfield (24250)
m$rMonExp_EGM[m$LGA_CODE11 == 26670] <- m$rMonExp_EGM[m$LGA_CODE11 == 24250]

#Moira(24900), Strathbogies(26430) and Gannawarra(22250) to Murrindini(25620)
m$rMonExp_EGM[m$LGA_CODE11 == 24900] <- m$rMonExp_EGM[m$LGA_CODE11 == 25620]
m$rMonExp_EGM[m$LGA_CODE11 == 26430] <- m$rMonExp_EGM[m$LGA_CODE11 == 25620]
m$rMonExp_EGM[m$LGA_CODE11 == 22250] <- m$rMonExp_EGM[m$LGA_CODE11 == 25620]

# merge fortifed spatial data with spatial data that now contains pokies data  
m2 <- merge(out, m, by.x="id", by.y="LGA_CODE11", all.x=TRUE)  
m2$exp <- cut(as.numeric(m2[,11]), breaks = c(seq(0, 12000, by = 2000)))  

# choose colours for the plot; we pick green for $$ 
library('RColorBrewer')
cbbPalette <- brewer.pal(6, "Greens")

library('ggplot2')
library('grid')

#ggplot code  to produce map 
map = ggplot(m2, aes(long, lat, group = group)) +  
  geom_polygon(aes(fill = exp), colour = I("black"), size = 0.2) +  
  scale_fill_manual(values=cbbPalette, name = "", labels = c("2k to 4k","4k to 6k", "6k to 8k","8k to 10K","10k to 12k")) +  
  ylab(NULL) +   
  xlab(NULL) + 
  theme(panel.background = element_blank(), plot.title= element_text(lineheight=2.8, face="bold"), 
        legend.text = element_text(colour="black", size = 12),
        legend.key.size = unit(.5, "cm"), 
        axis.line=element_blank(), axis.text=element_blank(), axis.ticks=element_blank(),
        plot.margin=unit(c(0,0,0,0),'cm'))
map

#output ggplot to .jpg 
g1 <- map + coord_map()
ggsave(plot=g1, filename=paste('pokies04to12V2', '.jpg', sep=''), scale=1.2,height=5, width=7)


